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We investigate the formation of a two-dimensional quasicrystal in a 
O ' monodisperse system, using molecular dynamics simulations of hard sphere 

' particles interacting via a two-dimensional square-well potential. We find 

that more than one stable crystalline phase can form for certain values of the 
g : square-well parameters. Qaenchmg the l.qa.d phase at a very low tempera- 

o 



ture, we obtain an amorphous phase. By heating this amorphous phase, we 
obtain a quasicrystalline structure with five-fold symmetry. From estimations 
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' of the Helmholtz potentials of the stable crystalline phases and of the qua- 
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■ sicrystal, we conclude that the observed quasicrystal phase can be the stable 
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phase in a specific range of temperatures. 



I. INTRODUCTION 

Stable quasicrystalline phases are typically found in binary mixtures [Q], where the vari- 
ous arrangements of the two components contribute to the degeneracy of the local environ- 
ments 0, allowing a quasicrystalline phase to be entropy stabilized p. With one notable 
exception previous studies did not support the existence of a stable quasicrystalline 
phase in a monodisperse system interacting with a simple potential 

We study a simple model that allows us to estimate the crystal and quasicrystal entropies 



and thereby study the Helmholtz potentials of the crystals and quasicrystal. The ground 
state of this system is a periodic crystal, yet we explore the possibility that the quasicrys- 
talline configuration is the equilibrium state in a certain temperature regime. Although 
quasicrystals do not have short range order, they do have recurring local environments that, 
in our model, resemble the basic cells of the stable crystalline phases. From the entropies 
of the stable crystalline phases and by estimating the configurational entropy of the qua- 
sicrystal, we infer that the quasicrystal may be an equilibrium state. We observe sharpening 
of five fold diffraction peaks when the starting amorphous phase is annealed. In two di- 
mensions, five fold diffraction peaks pertain to crystallographically disallowed point groups 
which characterize quasicrystals [0]. 

II. MD METHODS 

To study quasicrystalline stability in a monodisperse system, we perform molecular dy- 
namics (MD) simulations of a two-dimensional model of hard spheres interacting with an 
attractive square-well (SW) potential [Fig. Q]. The simplicity of this SW potential allows 
us to study the fundamental characteristics of the system. By tuning the width of the SW 
potential, we can control the local geometric configurations formed by the particles. The 
structures of the crystalline and quasicrystalline phases can thus be clearly defined and 
analyzed. 

We perform MD simulations in the NVT ensemble, using a standard collision event list 
algorithm |^ to evolve the system, while we use a method similar to the Berendsen method to 
achieve the desired temperature 0. The depth of the potential well is e = —1.0. Energies are 
measured in units of e, temperature is measured in units of energy divided by the Boltzmann 
constant, e/fcs, and the mass of the particle is m = 1. We choose the value of the hard core 
distance to be, a = 10, and the ratio of the attractive distance b to the hard core distance a, 
to heb/a = y/S. Since the diagonal distance between two corners of a square is \/2 times the 
length of one side, choosing b/a = \/3 favors the formation of a square crystal lattice where 
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each particle interacts with 8 neighbors [Fig. ^]. This constraint inhibits the formation of 
a triangular crystal, which would form at low temperatures if 6/a > y/S or at high densities. 

III. CRYSTAL AND AMORPHOUS PHASES 

Studying the behavior of the system at low temperatures we observe the formation of local 
structures similar to that shown in Fig. ^ These structures constitute local environments 
that can reproduce crystallographically allowed symmetry if translationally ordered. First, 
we consider the stable periodic crystal phases produced by translationally ordering each of 
the configurations in Fig. ^ and calculate the energies of these crystal structures at T = 0. 
In our system, the two allowed local configurations are the 4-particle square and the 5- 
particle pentagon (indicated by the symbol "P" in Fig. 0b, c). Particles form these two 
geometries because the nearest neighbor diagonal and adjacent distance between particles 
in these configurations is less than b/a = y/3, the SW width. Four particle squares make up 
the square crystal; since each particle has 8 neighbors, at T = 0, the potential energy per 
particle is Ugq = —4.0. Pentagons do not tile the plane; however, the formation of two kinds 
of crystals based on the local five particle pentagon is possible. In the type I pentagonal 
crystal, each crystalline cell consists of 5 particles, one of which has 8 neighbors and four of 
which have 9 neighbors; hence, Upi = — 4| [Fig. ||b]. In the type II pentagonal crystal, each 
crystalline cell consists of 6 particles, two of which have 8 neighbors and 4 of which have 9 
neighbors; hence, Upu = — 4| [Fig. |]c]. Since Upi < Upu < Usq, at our chosen density and 
low enough temperatures, the type I pentagonal crystal should be the stable phase at T = 
[Fig-i- 

Next, we investigate the stability of the three crystalline phases at T > by estimating 
the Helmholtz potential per particle A = U — TS in the square crystal and in the pentagonal 
crystals of type I and type II. Here S is the entropy. Since our simulations are performed at 
constant density, we must use the Helmholtz potential instead of the Gibbs potential. We 
study the system at dimensionless number density p = a^N/V = 0.857. We have simulated 
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a square crystal with N=961, a pentagonal crystal type I with N=1040, and a pentagonal 
crystal type II with N=792, all at the same p. We checked that at low temperatures, 
T < 0.1, the potential energy U{T) is temperature independent, and has the same value as 
the potential energy of the ideal crystal at T = 0. Hence, we approximate U{T) at higher 



In order to plot the behavior of the Helmholtz potentials of the three crystals for T < 0, 
we find the entropic contributions S", by estimating the entropy per particle for each of the 
three crystal types. We use the probability density y) to find a particle at position 
where the average is taken over every particle in the crystalline cell: 



We estimate y) by the fraction of the total time t spent by a particle in a discretized 
area, AxA?/, at a low enough temperature that the potential energy fluctuations of the 
crystalline structure are negligible. The values of the entropies for the three crystals are 
given in Table I. 

Our estimates for the temperature dependence of the Helmholtz potential for the three 
types of crystals are given in Fig. ^. The condition for stability of the pentagonal crystals 
is that their Helmholtz potentials, Apj and Apu, are lower than the Helmholtz potential 
of the square crystal, Asq. In accord with this condition, the square crystal is stable at 
temperatures above T = 0.203, the type II pentagonal crystal is stable between T = 0.195 
and T = 0.203 and the type I pentagonal crystal is stable below T = 0.195. 

While studying the interesting region around T ^ 0.2 (see Fig 0), we observe the for- 
mation of the quasicrystal. We choose to investigate, using MD simulations, our system 
at T ^ 0.2 because this is the temperature regime where the three crystals have similar 
values of Helmholtz potential. Cooling the fluid phase, we find the formation of the square 
crystal below T ^ 0.5. However, when further cooled into the temperature regime where the 
Helmholtz potentials of the two pentagonal crystals are lower than the Helmholtz potential 
of the square crystal, the system does not form pentagonal crystal I or pentagonal crystal II 
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(within our simulation times), but remains as the square crystal. Hence, we use a different 
approach to try to form the pentagonal crystals: we heat an amorphous phase. We first form 
the amorphous phase by quenching the system from high to very low temperatures T < 0.1. 
To do this, we study a system of = 961 particles at p = 0.857 which is initially 
in the fluid phase at high temperature, T = 10. We quench this system to T = 0.1 and 



thermalize for 10'' time units 1TT|. Time constraints prevent us from studying systems with 
more than 961 particles. Long thermalization times are required to stabilize thermodynamic 
observables like energy and pressure. 

The amorphous phase is a homogeneous mixture of pentagons and squares [Fig. The 
lack of long range structural order in the amorphous phase is evident from the homogeneity 



of the circles in the isointensity plot [Fig. When heating the amorphous phase 



to temperatures above T ^ 0.15, we find that diffusion becomes sufficient for local rear- 
rangement to occur, and the pentagons begin to coalesce. Instead of forming type I or II 
pentagonal crystals, the pentagons begin to form rows [Fig. ^] that bend at angles which 
are multiples of 36°. The angle in the bending of the rows gives rise to the five-fold orien- 
tational symmetry, which corresponds to the ten easily observed peaks in the isointensity 
plot [Fig. ^]. These 10 peaks are characteristic of the quasicrystal phase |]T3[, as they are 



arranged with disallowed fifth order point group symmetry 0. The configuration that we 
obtain has defects, mainly patches of square crystal, which cause the discontinuity in the 
rows and lead to the broadening of the diffraction peaks. For comparison, we present in 
Fig. 1^ the isointensity plots of the simulated square and pentagonal crystals. The diffraction 
patterns illustrate the symmetry of the original crystal system. The four equal sides of the 
square crystal unit cell [Fig ^a] are clear in the symmetry of the isointensity plot Fig. ^. In 
the isointensity plot of pentagonal crystal I [Fig. ^d], the central region which corresponds 
to the long range order, shows no hints of anything but well defined centered-rectangular 
symmetry [Fig ^] [jl4| . The isointensity plot of pentagonal crystal II has mainly a rectan- 



gular symmetry that matches the rectangular symmetry of the unit cells [Fig ^]. Although 
the two pentagonal crystals are formed from ordered pentagons, their long range symmetries 



are four sided. Their corresponding isointensity plots illustrate these four fold symmetries 
which are distinctly different from the five fold quasicrystal isointensity plot. 



IV. QUASICRYSTAL 



A. Formation 



Since the phase transition between the two pentagonal crystals occurs at T ^ 0.2, we 
choose this temperature as the one to investigate for quasicrystal formation. After the 
amorphous phase is quenched to T=0.1, we anneal the system at T = 0.205, for 2 x 10^ 



time units, and calculate the diffusion coefficient D, pressure P ||T^ and potential energy U. 



We calculate D using the Einstein relation D = ^lim^^oo ^ ^ where d is the system 

dimension. After a short initial period of increase, we observe that D and U decrease 
with time and reach plateaus [Fig. The diffusion coefficient approaches zero, which is 
consistent with the possible formation of a quasicrystal phase. The isointensity peaks also 
sharpen with the duration of annealing. Due to MD time constraints, we are not sure that 
we reach the potential energy of a perfect quasicrystal, which is expected to be comparable 
to the energies, Upi = — 4| and Upu = — 4|, of the pentagonal crystals. The lowest potential 
energy reached is Ugc = —4.25. 

We observe the spontaneous formation of the quasicrystal phase in the range of tem- 
peratures between T = 0.190 and T = 0.205. As we heat either the amorphous phase or 
the quasicrystal above T = 0.21 , the square crystal forms, consistent with the Helmholtz 
potential estimations of Fig. |. 

Next we address the question of whether the quasicrystal phase is stable, by comparing 
the values of the Helmholtz potential for the three crystal types. As can be seen [Fig. ^,d], 
the structure of the quasicrystal arises from the bending rows of pentagons which locally 
resemble the pentagonal crystals of either type I or II. We assume that local arrangements 



of particles corresponding to a square crystal are defects [16] that would be absent in the 
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perfect quasicrystal. If we assume that the local arrangement of the quasicrystal is similar 
to a combination of the local arrangements in the pentagonal crystal I and the pentagonal 
crystal II, we can approximate the Helmholtz potential of the quasicrystal by the average 
Helmholtz potential of the two pentagonal crystals. Because the quasicrystals have a positive 
entropy contribution to the total entropy due to their degeneracy we add an additional 
term —TSc to the original estimate of the Helmholtz potential energy. Here Sc is the entropy 
due to the possible configurations of the quasicrystal. 



B. Entropy 

We estimate Sc as the logarithm of the number of configurations formed by n pentagons 
in the quasicrystal. A single pentagon can be oriented in two possible ways when attached 
side by side to an existing row of pentagons. Neglecting the interaction between adjacent 
rows, we can estimate the upper bound for the number of configurations as 2", where n is the 
total number of pentagons in the quasicrystal. Note that at point A on Fig. ^, the Helmholtz 
potentials of both pentagonal crystals coincide, so an additional —TSc term should stabilize 
the quasicrystal in the vicinity of point A. 

To better estimate Sc, we notice that the bending rows of pentagons forming the qua- 
sicrystal resemble a compact self-avoiding random walk on the hexagonal lattice. The num- 



ber of such walks grows as where Z ^ 1.3 and n is the number of steps [0. Since 
the formation of one pentagon in the midst of a perfect square crystal lowers the energy 
of the system by f/ = — 1, we estimate n to be {Uqc — Usq)N. Assuming that the ground 
state energy of the quasicrystal is between Upj and Upu, the number of pentagons in the 
quasicrystal, should not be smaller than the number of pentagons in the crystal of type II 
(which is the pentagonal crystal with the lesser number of pentagons and has n = ^N). We 
estimate the entropy of configuration per particle to be Sc ~ \n{Z^)/N = | ln(1.3) = 0.087. 
Thus, the quasicrystal should be more stable than the pentagonal crystals between T = 0.16 
and T = 0.23 , where the gap between the Helmholtz potential of the pentagonal crystals 
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is smaller than the configuration term TSc which ranges from 0.014 to 0.020 in the interval 
where T increases from 0.16 to 0.23. Since the TSc term lowers the Helmholtz potential of 
the obtained quasicrystal configuration below the Helmholtz potentials of the two pentago- 
nal crystals, it is likely that the obtained state with five-fold rotational symmetry is not the 
coexistence of type I and II pentagonal crystals, but is a stable quasicrystalline phase. A 
more rigorous investigation of this problem would either require the construction of a perfect 



Penrose tiling [|1^,0] or of a random tiling [^,^ involving the local structures of crystals 
type I and II. 



V. DISCUSSION 

To summarize, perfect pentagonal crystals of type I and II do not form spontaneously 
during the time scales of our study. Instead, the quasicrystal, having long-range, five-fold 
orientational order with no translational order, forms from the coalescence of pentagons 
present in the starting amorphous phase. The starting amorphous configuration must ini- 
tially be quenched at a low enough temperature in order to prevent crystallization to the 
square phase. Moreover, the amorphous phase must be carefully thermalized at the quench 
temperature, as we have observed that, upon heating a poorly equilibrated amorphous phase 
with a higher concentration of squares, the system phase separates into regions of pentagons 
and squares. If the starting amorphous phase does not have a sufficient concentration of 
pentagons, the quasicrystal will not form: large regions of square crystal will inhibit the 
long range order of pentagons and thus not give rise to the 10 diffraction peaks in the isoin- 
tensity plot. It is interesting to notice that the bending rows observed in our quasicrystal 
could resemble the stripe structure of a spinodal decomposition 0. Anyhow, in the case of 
spinodal decomposition, the diffraction pattern would be similar to that of an amorphous 
structure. 

Before concluding, we note that Jagla, using Monte Carlo simulations, recently reported 
the existence of quasicrystals in a two-dimensional, monodisperse system of hard spheres 
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interacting with a purely repulsive potential [Q. The quasicrystal we observe has a different 
structure from that modeled by Jagla: our quasicrystal is not a ground state structure and 
forms only at nonzero temperature. Also, formation of quasicrystals in monodisperse systems 



has been observed using complex radially symmetric potentials both in two dimensions ||20 
and three dimensions P^ , |25|] . To the best of our knowledge, the quasicrystal found in our 
simulations has a structure different from those previously studied. 
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TABLE I. Energy U, entropy S and the Helmholtz potential A at temperature T = 0.2 where 



the quasicrystal is found 



Crystal 


U 


S 


A{T = 0.2) 


Pentagonal I 




1.259 ±0.028 


-4.652 


Pentagonal II 


-4 


1.603 ± 0.0052 


-4.654 


Square 


-4 


3.247 ±0.021 


-4.649 
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FIG. 1. The square- well potential. The ratio of the attractive distance b to the hard core 
repulsive distance a is b/a = \/3. The depth of the square-well e = —1.0 is the interaction energy 
per pair of particles. 
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FIG. 2. Repeating segments of the three crystals, (a) In the square crystal, each particle interacts 
with 8 nearest neighbors. The square crystal is constituted by particles interacting with a square 
geometry, (b) In type I pentagonal crystals, ^ of the particles have 8 neighbors and | of the 
particles have 9 neighbors, (c) In type II pentagonal crystals, | of the particles have 8 neighbors 
and I of the particles have 9 neighbors. Five particle pentagons, denoted by letter "P" , form the 
pentagonal crystals. The particles indicated in white are the particles in a basic cell that can be 
used to construct the crystal by translation; there are respectively one, five and six particles in the 
unit cell of the square, pentagonal I and pentagonal II crystals. 
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FIG. 3. The Helmholtz potentials of pentagonal crystals of type I and II and the square crystal at 
various temperatures. Points A, B, and C, of the inset, indicate the intersections of the Helmholtz 
potential lines at Ta = 0.195 ± 0.010, Tb = 0.201 ± 0.005, Tc = 0.203 ± 0.006. The solid line 
indicates the lowest Helmholtz potential: below Ta the type I pentagonal crystal is the most 
stable, between Ta and Tc the type II pentagonal crystal is the most stable, and above Tc the 
square crystal is the most stable. 
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FIG. 4. Amorphous and quasicrystal phases are shown along with their corresponding isointen- 
sity plots: the simulated equivalent to a crystallographic diffraction pattern, given by the Fourier 
transform of the density function: the darkness is proportional to the amplitude of the Fourier 
transform. On the original system snapshots (a,c), pentagons are indicated by the shaded areas 
and lines indicate interacting pairs of particles, (a) Uniformly distributed pentagons in the amor- 
phous phase give rise to the (b) homogeneous rings in the isointensity plot, (c) The pentagons in the 
quasicrystal phase have coalesced in curved rows that run approximately parallel to one another, 
in contrast to part (a) where the rows are much less apparent and are not even approximately 
parallel, (d) The ten isointensity peaks of the quasicrystal. 
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FIG. 5. The (a) square crystal, (b) type I pentagonal crystal and (c) type II pentagonal crystal 
isointensity plots. 
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FIG. 6. Behavior of pressure, P, potential energy per particle, U and diffusion coefficient, D 

versus time when the system, initially in the amorphous phase, is equilibrated at T = 0.205. The 
density is p = 0.857 and the number of particles is A'" = 961 
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